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Abstract 

A description in terms of phase and amplitude variables is given, for nonlinear oscillators subject to white 
Gaussian noise described by Ito stochastic differential equations. The stochastic differential equations derived 
for the amplitude and the phase are rigorous, and their validity is not limited to the weak noise limit. If 
the noise intensity is small, the equations can be efficiently solved using asymptotic expansions. Formulas 
for the expected angular frequency, expected oscillation amplitude and amplitude variance are derived using 
Ito calculus. 
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1. Introduction 


Nonlinear oscillations are ubiquitous in natural sciences and technology [11,0. An ideal oscillator 
would exhibit a perfectly periodic behavior, represented by a limit cycle in the state space. However, the 
output of actual oscillators is always corrupted by different types of disturbances, such as internal noise 
sources, thermal noise, interactions with the environment and with other systems. The autonomous nature 
of oscillators implies that any time shifted version of a solution is also a solution, therefore at least one 
direction exists along which perturbations are preserved. This explains why oscillators are very sensitive 
to noise, and introducing even small random perturbations into oscillators leads to dramatic changes in 
their frequency spectra and timing properties. This phenomenon, peculiar to oscillators, is known as phase 
noise or timing jitter. As a consequence, characterizing how noise affects the dynamics of oscillators is of 
paramount importance. 

Phase models and phase reduction methods are powerful tools for analyzing the effect of perturbations 
on oscillators 00 , 0 , 0 - Phase models are based on the idea to describe the state of an oscillator using a 
phase and an amplitude functions. The phase represents the projection of the perturbed trajectory onto a 
reference orbit, e.g. an unperturbed limit cycle. The amplitude function, instead, measures the deviation of 
the perturbed orbit from the reference limit cycle. This picture can be further simplified under the hypothesis 
that the the limit cycle is asymptotically stable, and the noise is weak. With these assumptions, one expects 
that deviations from the limit cycle remain small, and that the amplitude function can be approximated by 
its unperturbed value. Then only the phase dynamics is retained, obtaining a phase reduced model. 

In the last few years, phase models and phase reduced models have been extensively studied for the 
systematic investigation of the phase noise problem in oscillators MiSEJEl- On the one hand, phase 
reduction method proved to be a valuable tool to derive simplified, mathematically tractable models, whose 
analysis has significantly increased our understanding of the influence of noise on oscillators. On the other 
hand, most of these models are of little help in analyzing practical oscillators. In fact, an explicit formula 
for the phase function can be found only for few, trivial, oscillators. In all other cases one must resort to 
numerical methods that are either approximate in nature, or unsuitable for oscillators of order higher than 
the second 12|, 13, 14, 10. The consequence is twofold. First, phase models cannot be used to study and 
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design practical oscillators, and second, experimental and/or numerical results obtained from real world 
oscillators cannot be readily used to test the validity of phase models. Another potential problem concerns 
the model order reduction procedure. Amplitude fluctuations can be neglected only if the relaxation to the 
stable orbit is instantaneous, or at least, if it occurs on a time scale much shorter than the typical time 
for the phase dynamics. Such an assumption is often taken more for mathematical convenience than being 
physically plausible. Moreover, when dealing with noisy perturbations the stochastic nature of the processes 
should be taken into account. Appropriate correction term should be introduced before__order reduction, to 
deal with the possible correlation between stochastic variables and noise increments 


lflllZl- 


This paper proposes possible solutions for the aforementioned problems. A novel amplitude and phase 
description for nonlinear oscillators subject to white Gaussian noise is presented. The method is based on 
the generalization of a classical technique 18| to the case of Ito stochastic differential equations. Contrary 
to other derivations, the proposed amplitude and phase description does not rely on undetermined phase 
functions, and its validity is not necessarily limited to the case of weak noise. It is shown that using an 
appropriate basis, a partial decoupling between the amplitude and the phase dynamics can be obtained. The 
resulting equations represent the ideal starting point for the derivation of phase reduced models. It is also 
shown that under suitable simplifying hypothesis, the proposed model reduces to other models previously 
described in literature, i.e. previously proposed models are special cases of our description, obtained applying 
different degrees of approximation. For the case of weak noise, the amplitude and phase equations can be 
analyzed using an asymptotic expansion method. The equations are recast as a sequence of linear stochastic 
differential equations, that can be solved iteratively. Instead to solve these equations, we exploit the linearity 
and the properties of Ito integrals to obtain a statistical characterization for the noisy oscillators, finding 
analytical formulas for the expectation values of the angular frequency, the amplitude and the amplitude 
variance. As an example, the technique is applied to a noisy Stuart-Landau oscillator. 


2. Amplitude and phase dynamics of noisy oscillators 


Nonlinear oscillators subject to white Gaussian noise can be conveniently described by the stochastic 
differential equation (SDE) 

dX t = a(X t ) dt + eB(X t ) dW t (1) 

where X t : 1R, H > R ra is a stochastic process describing the state of the oscillator, a : R n i-> R™ is a vector 
valued function that defines the oscillator dynamics, B : R ra R" ,m is a real valued n x m matrix, e is a 
parameter that measures the noise intensity and W t : R i —> R m is a vector of Brownian motion components. 
We shall assume a £ C k (fl. C R"), with k > 1 and that B satisfy a Lipschitz condition, to guarantee the 
existence and uniqueness of the solution M- 

Depending on the definition adopted for the stochastic integral, the SDE |[lj) can be interpreted following 
two main schemes: Stratonovich or Ito. If Stratonovich interpretation is used, the amplitude and the phase 
equations can be derived in a straightforward way using the standard procedure described in [l^|, since 
traditional calculus rules apply. However the analysis of the resulting equations gets more difficult because 
of the “look in the future property” of Stratonovich stochastic integral [19| ■ By contrast Ito interpretation 
does not suffer of the anticipating nature, but a new set of calculus rules, known as Ito calculus, must be 
used. Therefore, we shall assume that CD is an Ito equation, and we shall use Ito calculus to derive and to 
analyze the amplitude and phase equations. 

For e = 0 equation CD reduces to the system of ordinary differential equations (ODEs) 

A/) = «<*(()) (2) 

We assume that © admits a T-periodic solution, i.e. a function x s (t) exists that satisfies © with the 
property x s (t + T) = x s (t). Before proceed with the main result we introduce some notation. We define the 
unit vector tangent to the limit cycle (the symbol ' denotes the derivative with respect to the argument) 


Ul(t) = 


x s(t) 

\ X s(t)\ 


( 3 ) 
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Together with ui(t) we consider other n — 1 linear independent vectors 1 x 2 (t ),..., u n (t ), such that the set 
{■Ui(i),... ,u n (t)} is a basis for R™, for all t. We remark that, differently from traditional derivations [18| , 
we require the vectors {u\(t ),..., u n (t)} to be linearly independent, but not necessarily orthogonal. The 
conditions for the existence and a procedure to construct an orthogonal set are discussed, for instance, in 
00- Together with {ui(t), ..., u n (t)} we also consider another basis, {iq (t),... ,v n (t)} constructed 
as follow: Given the matrix U(t) = [tti(t),..., we define the reciprocal vectors vf (t ),..., vf(t) 

to be the rows of the inverse matrix V(t) = t7 _1 (t). Thus {iq(t),... ,v n (t)} also span R ra and the bi 
orthogonality condition vf Uj = uf Vj = 5^ holds. We shall also use the matrices Y(t) = [u 2 (t ),..., u n (t)], 
Z(t ) = [v 2 (f),... ,v n (t)\, and the modulus of the vector field evaluated on the limit cycle, r(t) = |a(® s (t))|. 

The most important concept to be defined in the analysis of oscillator noise is the phase concept. A 
phase function is intended to represent the projection of the oscillator’s state onto a reference trajectory, 
normally the unperturbed limit cycle. Being associated to a neutrally stable direction, random fluctuations 
of the phase will be neither adsorbed nor amplified, they persist and may eventually accumulate with time 
growing unboundedly large. This explains why the phase is a privileged variable. 

We introduce a phase function 6 : R™ 1 —> [0,T), interpreted as an elapsed time from an initial reference 
point. Consider a point aq(0) on the limit cycle, and assign phase zero to this point, i.e. #(a; s (0)) = 0. The 
phase of the point x s {t) is 6(x s (t)) = t, mod T. Thus, the phase represents a new parametrization of the 
limit cycle. Together with the phase function we shall consider an amplitude function R : R n R n_1 , with 
9,R £ C m (R n ), m > 2. The amplitudal function R(x) is interpreted as an orbital deviation from the limit 
cycle. The following theorem represents the generalization of a classical result [18] , to the case of Ito SDEs. 

Theorem 1. Consider the ltd SDEs © such that the ODEs obtained setting e = 0 admit a T-periodic 
limit cycle x s (t). Let {ui(t),... ,u n (t)} and {v\(t),... ,v n (t)} be two reciprocal bases such that U\{t) sat¬ 
isfies © and such that the bi-orthogonality condition vj Uj = uf Vj = Sij holds. Consider the coordinate 
transformation 



Then a neighborhood of the limit cycle x s (t) exists, where the phase 9{t) and the amplitude R(t) are ltd 
processes and satisfy 


d9=[ 1 + a\ ( 9 , R) + e 2 ai(0, R)] dt + eB x (9, R) dW t 
clR=[L{9)R + a 2 {9 , R) + e 2 a 2 (6», R)\dt + eB 2 (9 , R) dW t 


( 5 ) 


( 6 ) 


1 We shall use the term “amplitude” instead of the more correct “amplitude deviation” for the sake of simplicity. 
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with (explicit dependence on 9 and t is omitted for simplicity) 
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(9) 
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( 11 ) 

( 12 ) 

(13) 


Proof: First we show that a neighborhood of the limit cycle exists, where 9 and R are Ito processes. To 
this end consider the coordinate transformation (U) and the associated Jacobian matrix 


Dh(9,R) = 

On the limit cycle R = 0 and then 


8h 8h 

~89 8R 


= [x' s {9) + Y\9)R Y{ 


(14) 


Dh(9M)\ R=0 = [x' s (9) Y(9)} = (r(9) Ul (9),u 2 (9),...,u n (9)} 

Since {tti(f),... ,u n (t)} is a basis for ]R n , it follows that the determinant of the Jacobian matrix is not 
null. Then by the Inverse Function Theorem there exists a neighborhood of R = 0 where h is invertible. 
Moreover, if h is of class C k then its inverse is also of class C k . Taking the inverse of h we can write 9 = 9(x) 
and R = R(x ), and if the basis vectors are smooth enough it follows from Ito formula that 9 and R are Ito 
processes. 

Using Ito formula and eq. m, x = ft.(0, R) implies 


8h 


8h 


1 8 2 h 


1 


8 2 h 


dx d9 + f^dR + - ^i(d9) 2 + r 2 dR T ^ 
=a(h{0, R))dt + £ B(h(6, R)) dW t 


dR T — 


1 8 2 h 


2 898R 


d9dR 


(15) 


where 8h/8R and 8 2 h/8R 2 are the matrices of first and second partial derivatives with respect to the 
components of R , respectively. Introducing (Hj) in (fl5l) yields 


8Y 


89 


a(x s ) + — R)d9 + YdR+ - 


1 / 8a{ x s ) 8 2 Y 


2 V 89 


89 2 


8Y 

R)(d9) 2 + — d9dR 


= a(x s + Y R)dt + eB{x s + YR)dW t 


(16) 


Multiplying to the left by vj and using the bi-orthogonality condition we get 


,8Y 

~89 


R\d9 


1 


8a(x s ) 8 2 Y 


89 


+ -Qp R )WY+vi 


8Y 

~89 


d9 dR 


= vja(x s + YR)dt + £Vi B(x s + YR)dW t 


(17) 
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Figure 1: Two possible decomposition of the stochastic process Xt- At the time t\ the process is decomposed as Xt 1 = 
x s (0(ti)) + Y(Q(ti))R(ti) using two different basis vectors. Left: orthogonal basis. Right: “oblique” basis. Red line is the 
stochastic process, blue line is the limit cycle shown for reference. 


By converse, multiplying ci to the left by Z T yields 


.dY 


Z 1 -RdO + dR + -Z 


1 T f 9a ( x s) . d 2 Y 


89 


,d Y 


+ ^ R )^ 2 + zT ie dddR 


2 V 80 

= Z T a{x s + YR)dt + eZ T B(x s + YR)dW 


(18) 


Since 0 and R are Ito processes they satisfy relations of type dO = a dt + e/3 dW t , and dR = 7 dt + ecr dW t , 
respectively. By Ito lemma ( dO ) 2 = £ 2 f3f3 T dt, dO dR = £ 2 ct(3 t dt. Introducing these results in (fTTl) . (1151) 
and equating terms in dWt we obtain 


P = (r + vJ^Rj vf B(x s + YR) 

dY 

cr =Z T B(x s + Y R) — Z T R/3 

Ou 


(19) 


( 20 ) 


Finally, using (fTHD . (d, together with ( dO ) 2 = (3(3 T dt, dO dR = a(3 T dt in (fTTl) . d and rearranging the 
terms we get the thesis. □ 

Figured shows the idea behind the decomposition defined by (d- At any time instant ii, the stochastic 
variable X tl is split into two components. One component is tangent to the cycle, and it is given by the 
limit cycle evaluated at the stochastic time 0(ti). The second component is transversal to the cycle, it is 
represented by the distance R(ti) measured along the vectors U2(0(ti)),..., u n {0(t.\)). Although the choice 
of an orthogonal frame may look the most natural, the use of particular non orthogonal frames offers some 
advantages, as we shall show in the next section. 


3. Amplitude fluctuations and phase—amplitude decoupling 

The amplitude and phase equations m and © depend upon the choice of the basis vectors {tti(f), ... ,u n {t)}. 
It is natural to ask whether a basis exists that has to be preferred to others. In particular, we shall show 
that, if the basis {ui(t),..., u n {t)} is chosen according to Floquet’s theory, the phase dynamics can be par¬ 
tially decoupled from the amplitude dynamics. First we show that, under suitable conditions, the amplitude 
variable almost surely remains in a neighborhood of the unperturbed limit cycle. The following theorem 
represents an adaptation of a classical result [ 11 ], to the case of bi-orthogonal basis. 


Theorem 2. Consider the noiseless oscillator m with the T-periodic limit cycle x s (t). Let 1,/X2, 
he the characteristic multipliers of the variational equation 


dy(t) 

dt 


= A{t)y(t ) 

5 


> A In 


( 21 ) 





where Aft) = da ^^ t '>) j s Jacobian matrix evaluated over x s (t). Then R = 0 is an equilibrium point for 
the amplitude equation with characteristic multipliers /x 2 , • • •, y n - 


Proof: For e = 0 the amplitude and phase equations ©, ©, reduce to 

dd 


dt 

dR 

dt 


—1 + a\(9 , R) 
=L{d)R + a 2 {9,R ) 


( 22 ) 


(23) 


where ai(d, R), a 2 (d,R) and Lid) are given by (0, (flOl) and m, respectively. It is trivial to verify that 
ai(9, 0) = 0 and a 2 (9, 0) = 0, thus R = 0 is an equilibrium point for the amplitude equation (H?l) . The 
variational equation for the amplitude is obtained by taking the Taylor expansion of a(x s + YR) around 
x s (9) in the equation for a 2 (0, R) (eq. (fill) I. Using the fact that Z T a(x s ) = 0 and neglecting 0(R 2 ) terms 
we get 


dR 

dt 


L{9) + Z T AY 


R 


(24) 


Next we show that if R s (t) is a solution of (HO) , then y s (t) = ax' s (9) + Y ( 9)R s (t ) with a £ R, is a solution 
of (l2ll) . The following equality must hold 


dy s = d 2 x s 
dt a 89 2 


dY 

~d9 


R s 


, dR s 

dt 


= A a 


dx s 

dd 


YR, 


(25) 


where (1221) has been used and higher order terms in R have been neglected. It is well known that if x s 
solves ©, then x' s solves (ED- Rearranging the terms, equation (HZHll reduces to 


Y dR 1 = _dY R 

dt dd s 


AYR , 


Multiplying to the left for Z r , and using the bi-orthogonality condition we obtain 


(26) 


(27) 


that coincides with (1241) . Let Q(t) be a fundamental matrix solution of the amplitude variational equation 
EH, then a fundamental matrix solution of the variational equation ED is of the form 


*(<) = 




(28) 


where [a ±,..., a n —i\ is a row vector of real constants. Since (1211) is a linear system with T-periodic coeffi¬ 
cients, there exists a constant matrix C with eigenvalues 1 ,^ 2 , ■ ■ •, £t n , such that 4>(f + T) = 4>(t)C. From 
the periodicity of x s it follows that C has the structure 


C = 


1 c 2 

0 Cl 


(29) 


where the eigenvalues of Ci are y 2 ,...,y n . On the other hand <h(f + T) = $(t)C together with the 
periodicity of d imply 


OX ux 

qq [cx-i, ■ ■ ■, ot n - 1 ] + YQ(t + T) = -QQ - (C 2 + [ai,..., a n ]Ci) + YQ(t)C i 


(30) 


Multiplying to the left by Z T we obtain Q(t + T) = Q(t)C±. Then the fundamental matrix solution of the 
amplitude variational equation (l24l) has multipliers fi 2 ,..., y n as required. □ 
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Theorem [2] implies that the amplitude equilibrium point R = 0 inherits the stability property of the limit 
cycle x s (t), e.g. if x s (t) is asymptotically stable (i.e. /X 2 , • • ■, Abi have modulus less than one), then R = 0 
is also asymptotically stable. If this is the case, we expect that trajectories of the perturbed trajectory 
remains in a small neighborhood of the unperturbed limit cycle. 

Since the amplitude is robust to perturbations, in most practical applications the phase is the only 
relevant variable. The following theorem establishes that using Floquet’s theory, a particular basis can 
be identified such that a partial decoupling of the phase dynamics from the amplitude dynamics can be 
achieved. Before introducing the theorem we recall that from Floquet’s theory, the fundamental matrix 
solution of the linear time periodic variational equation m can be written in the form 

m = P(t)e Dt S 0 (31) 

where P(t) is a T-periodic matrix, So = P - 1 (0), and D = diag[iq,..., i/„] is a diagonal matrix whose 
diagonal entries are the Floquet’s characteristic exponents Ej |. 

Theorem 3. If the basis vectors {u\(t ),... ,u n (t)} are chosen such that [r(t)ui(t), u 2 (t ),..., u„(t)] = P(t) 
then, the ltd processes for the phase and amplitude reduce to 

dd =(1 + oi(0, R) + £ 2 di(0, R))dt + eB^O, R) dW t (32) 

dR =(D R + 0 , 2 ( 0 , R) + s 2 0 , 2 ( 0 , R))dt + eB 2 (6, R) dW t (33) 

where D = diag[v 2 , • ■ ■, Vn\, and the Taylor series ofdi(0, R), a 2 (0, R), do not contain linear terms in R. 

Proof: It is well known from Floquet’s theory that the columns of P(t) are linearly independent for 
any t, and thus they can be chosen as a basis for R ra . Moreover, x' s (t) is a solution for the variational 
equation, associated with the structural Floquet’s exponent vi = 0. Thus we can take the first column of 
P(t) to be r(t)ui(t) = x' s (t). Equation (13T1) implies P(t) = <&(f)P(0) e~ Dt . Taking the derivative of (l3TT) 
we have $'(i) = P' e Dt So + P D e Dt So, and taking into account that d>(t) is a fundamental matrix of the 
variational equation (l2TTl yields P' = AP — PD. Removing the first column we get Y' = AY — YD. 
Substituting for Y' in 0 ,' m and ED. taking the Taylor series o(x s T -R( — o(x s ^ T AA^R ..., and 
using the bi-orthogonality condition the thesis follows. □ 

The partial decoupling established by theorem [3] together with the fact that, under the asymptotic 
stability hypothesis the amplitude fluctuations remain small, permit to derive simplified phase models. In 
particular, if amplitude fluctuations are very small, an adiabatic approximation R(t) = 0 may be used. If 
0(e 2 ) terms are neglected, the phase dynamics reduces to 

d0 = dt + £B 1 (0,O)dW t (34) 

that is equivalent to the traditional phase model derived in 0 , 0 . Conversely, retaining 0(s 2 ) terms yields 
the phase equation 

d0= (l+£ 2 a 1 (0,0))dt + EB 1 (0,0)dW t (35) 

analogous to those derived in jl 6 |, [TVi]. 


4. Small noise limit and asymptotic expansion 

In general the phase and amplitude equations 0 and® are not easier to solve than the original SDE ®. 
On the other hand, the phase reduced model (15H) leads to incorrect predictions, such as that the noise does 
not influence the expected angular frequency of the oscillator 16]. Finally although it is more consistent from 
the point of view of physics, the reduced phase equation © may lead to inaccurate results as a consequence 
of the adiabatic approximation. In fact due to the nonlinear nature of oscillators, perturbations along certain 
directions are amplified, while other are reduced. The result is a net contribution to the amplitude deviation 
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so that, contrary to the adiabatic assumption, the expected value of the amplitude deviation E[R] is not 
null. 

In what follows an alternative approach is proposed, which is based on the simultaneous solution of 
the phase and the amplitude equations. In the weak noise limit (e -C 1) equations ([5]), and ©, can be 
efficiently solved using asymptotic expansions. For the sake of simplicity, we restrict the attention to second 
order oscillators, so that the amplitude deviation is a scalar variable. Higher order oscillators do not pose 
any particular problem, they only make the notation more involved. We search for solutions in the form 

6 = 9q + edi + e 2 9 2 + ■ .., R = Rq + eR\ + e 2 R 2 +_ Introducing these ansatzs in (J5]) and (0 , and equating 

the same powers of e we obtain: 
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89 
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(41) 


In equations (1351) (TUI) . L and its derivatives are evaluated at 6*o, while ai, a 2 , ai, a 2 , B i, B- 2 , and their 
derivatives are evaluated at (9q,Rq). The initial conditions are conveniently chosen in the form 9 n ( 0) = 0, 
ii n (0) = for all n. The zeroth order equations (l36l) . (l37l) are nonlinear ODEs describing the noiseless 
oscillator. They have the simple solution 9q = t, Rq = 0, that represents the limit cycle x s (t) in the 
phase and amplitude variables. Starting with equations (1551) (1551) . we have a sequence of time dependent 
Ornstein Uhlenbeck processes, whose analytical solution can be found using the method of integrating 
factors 0- The first order solutions are introduced into second order equations, which can be solved by the 
same method. The process can be iterated to higher powers of £ to obtain increasingly accurate solutions. 
The power series solution will not normally be a convergent series. However it is possible to show that 
the expansion is asymptotic, in the sense that the difference between the exact solution and its e n order 
approximation is of order e n+1 [ri . 

The explicit solution of equations (1551) (Till) is of little use, because it depends on the particular realization 
of the stochastic process Wt, he. different realizations of the Brownian motion lead to different specific 
solutions. Nevertheless useful information can be obtained without solving the SDEs. Taking the expectation 






















values on both sides of (1381) (HT1) and using the zero expectation property of Ito integral, the following ODEs 
for the expectation values are found 


dE{91 
dt 

dE[R i 
dt 
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~0R 
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(42) 

(43) 
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dt 
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d 2 a\ 

' 090R 
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(44) 


(45) 


Due to the linearity of the SDEs (1551) (HU) , system (H51) (Tl5l) can be closed. Using Ito formula, SDEs 
for d{9 2 ), d{R\) and d(9\ R\) can be derived, and taking stochastic expectations the following additional 
equations are found 


dE[9 2 } 

dt 

dE[R\] 

dt 

dE[9 1 i?i] 
dt 


=2 ^E[9\] + ^E^R,} + B.Bl 

=2(i + ^) E i R i] + 2 (|^° + ^f) E ^ + B ^ 

< L - + w + m w ^ + m E ^ + ^ + + B ' B * 


(46) 

(47) 

(48) 


Equations tEj m form a nonhomogeneous linear system of ordinary differential equations with time 
periodic coefficients. The equilibrium points allow to calculate the stationary expected angular frequency 
E[d9/dt] = 1 + sE[d9i/dt] + £ 2 E[d9 2 /dt], the expected amplitude E[R] = eE[Ri] +e 2 E[R 2 ) and amplitude 
variance var(i?) = e 2 var(i?i) + 0(e 3 ). 


5. Application 


As an example of application, we consider a Stuart-Landau oscillator with multiplicative noise 

dip = (a — (3p 2 ) dt + e pdWi 
dp = (p — p 3 ) dt + £ p 2 dW 2 


(49) 


In absence of noise (for £ = 0), the Stuart- Landau system admits an asymptotically stable limit cycle 


x s {t) 


(a - 0)t 

1 


(50) 


so that the unit tangent vector is U\(t) = [1,0] T . 
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5.1. Phase and amplitude equations with orthogonal basis 

Consider the orthogonal basis composed by ui(t) = [1,0] T , and u 2 {t) = [0, — 1] T . Obviously Vi(t) = tij(t) 
for i = 1,2. With this basis, the change of coordinates x(t) = h(6 1 R) = x s {9) + u 2 {9)R implies <p = ( a — p)9 
and p = 1 — R. Using equations 0 m it is straightforward to derive the phase and amplitude equations 


de = i l + r^( R - R2 )) dt + T^\( l - R ) dW i 
V \ a — P\ ) \a-p\ 

dR = (—2i? + 3i ? 2 — R 3 ) dt — £ (1 — R) 2 dW 2 


As expected in the equation for the phase a drift term linear in R appears. The equations for the first few 
terms of the asymptotic expansion are 


dd o 



2 /? 

W~P\ 


(■Ro - R 2 0 ) 


dt 


(52) 


dRo = (—2i ? 0 + 3i?o — Rq) dt 


(53) 


2 3 1 

dd i =1 —^ (1 - 2Ro) i?i dt + T -- [(1 - R 0 ) dW x (54) 

F-P| |a-/3| L 

dRi = (-2 + 6R 0 - 3 R$) R 1 dt - (1 - R 0 ) 2 dW 2 (55) 

23 1 

de 2 [ (1 - 2i? 0 ) R 2 - Rf\ dt - -- R 1 dW 1 (56) 

|a-/3| |a-P| 


dR 2 = [( 6 i ? 0 - 3i?p - 2) #2 + 3(1 - R 0 )R 2 1 \ dt + 2(1 - R 0 ) R 1 dW 2 (57) 

Taking the stochastic expectation on both sides of (l54l) ([57]) and using the zero expectation property of Ito 
stochastic integrals yields 


dE[9\] _ E ' dO i 
dt dt 


—(1 — 2Ro)E[Ri] 
ol-P | 


(58) 


dt 

dE[9 2 } 
dt 


= (-2 + 6R 0 



- 3R 2 ) ElRj] 

^l^[-E[R 2 1 ] + (l-2R 0 )E[R 2 ]] 


(59) 

(60) 


= ( 6 R 0 - 3R 2 0 - 2) E[R 2 ] +3(1- Ro) E[R\] (61) 

To close system (1551) (lOll) . an ODE for E[Rf\ is needing. Using Ito formula the SDE d(i?i ) 2 = 2RidR\ + 
(dRi) 2 is obtained, with ( dR\) 2 = (1 — Ro) 4 dt as a consequence of (l55l) and of Ito lemma. After substitution 
and taking the expectation value the closing equation is found 

= (-4 + 12i ? 0 - 6 i?g) E[R\] + (1 - i ? 0 ) 4 (62) 

Integrating m by separation of variables it is found that, as expected, Ro{t) —► 0 asymptotically for 
t —> +oo. Taking this into account the asymptotic stationary solution of system (I5%1) (Ifnil) are easily found 
E[Ri] = 0, E[R\] = 1/4, E[R 2 ] = 3/8, E[d9\/dt\ = 0, and E[d9 2 /dt\ = /3/(4|a — P |). It follows that the 
expected angular frequency, amplitude and amplitude variance are 


E 


'd9' 

dt 


1 + c 


p 

4|<a — P\ ’ 



Var [R] = 



(63) 
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5.2. Phase and amplitude equations with Floquet’s basis 
The Jacobian matrix evaluated over the limit cycle is 


A(t) 


0 -2/3 

0 -2 


(64) 


with eigenvalues pi = 0, p 2 = —2. The associated eigenvectors are the Floquet’s vectors u\(t) = [1,0] T 
and u 2 (t) = [/3, 1] T , while inverting the matrix U{t) = [ui(t), u 2 (t)\ we find the Floquet’s co-vectors 
vi(t) = [1,-/3] and v 2 {t ) = [0,1]. Repeating the calculations of the previous section we find that the 
relation between the old and the new coordinates is ip = (a — (3)6 + f3R, p = 1 + R. The phase and amplitude 
equations in the new basis are 

d9 = { 1 + ]^\ {2 + R ) ) dt + ]^ W\ {1 ~ R)dWl ( 65 ) 

dR = - (2 R + 3 R 2 + R 3 ) dt + e{ 1 + R) 2 dW 2 

Conversely to the orthogonal case, the phase equation obtained using Floquet’s basis does not contain a drift 
term linear in R, in accordance with theorem [3] The new equations for the first terms of the asymptotic 
expansion are 


dd o — 


1 + 


P 

|a-/3| 


(2 Rl + Rl) 


dt 


( 66 ) 


dR 0 — — (2jRq “I - 3i?Q + R 3 0 ) dt 


(67) 


dd 1 (4i?o + 3 Rl) R x dt + [(1 + Ro) dW x - /3(1 + R 0 ) 2 dW 2 ] (68) 

\a — p\ \a — p\ L 

dRi = - (2 + 6i? 0 + 3i?o) Ri dt + { 1 + R 0 ) 2 dW 2 (69) 

dd 2 =r dL— [(4i? 0 + 3 R 2 )R 2 + (2 + 3R 0 )R 2 1 ] dt 

\a — p\ 

+ [-Rl dW i ~ 2/3(1. + Ro)RidW 2 ] (70) 

|a-/3| L 


dR 2 = - [(2 + 6R 0 + 3 Rl)R 2 + 3(1 + i? 0 )i?i] dt + 2(1 + R 0 )R 1 dW 2 
The equations for the stochastic expectations read 


dE\di] 

dt 

dE[Ri\ 

dt 

dE[d 2 \ 

dt 

dE[R 2 } 
dt 

dE[R\\ 


=E 


ddi 

dt 


P 


r (4i2 0 + 3iZg)f;[i2 1 


= - (2 + 6i? 0 + 3 R 2 0 ) £?[i2i] 


=E 


dth 

dt 


P 


■ [(4i? 0 + 3 Rl)E[R 2 ] + 2 E[R - l 2 ] + 3 R 0 E[Rl}] 


= - (2 + 6i? 0 + 3i?g) E[R 2 ] - 3(1 + Rq)E[R 2 ] 
= - (4 + 12i? 0 + 6i?g) E[R 2 } + (1 - R 0 ) 4 


(71) 

(72) 

(73) 

(74) 

(75) 

(76) 
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Figure 2: Accuracy of the asymptotic expansion: phase and amplitude error vs time for various order od approximation and a 
specific realization of the Brownian motion. Floquet’s basis is used. Parameters are a = 5, = 1, and e = 0.15 



Figure 3: Left: Expected normalized angular frequency vs noise intensity, for the Stuart-Landau oscillator. Right: Expected 
amplitude deviation vs noise intensity. Theoretical predictions are obtained using both the orthogonal basis and Floquet’s 
basis. The values of the parameters are a. = 5, /3 = 1. 


where equation (EH) has been obtained by the same procedure described in the previous section. The 
asymptotic stationary solutions are E[Ri] = 0, E[R\] = 1/4, E[R 2 ] = —3/8, E[d9i/dt\ = 0, E[d9 2 /dt\ = 
f3/(2\a — /3|). Finally, the expected angular frequency, amplitude and amplitude variance are expected 
angular frequency are 


E 


'd9' 

dt 


1 + e 2 


P 

2|<a — f3\ ’ 



Var[l?] = -e 2 


(77) 


The theoretical predictions have been compared with Monte-Carlo simulations. The Stuart Landau 
equation (1491) has been integrated numerically using both Euler-Maruyama and Milstein integration schemes. 
Figure [2] shows the difference versus time, between the numerical solution of (HiH) and the asymptotic 
expansion solution, for different order of approximation. The result has been obtained using Floquet’s basis 
and is relative to a specific realization of the noise. It is indicative of the increasing level of accuracy of the 
expansion. 

Figure [3] shows the expected amplitude and angular frequency versus the noise intensity, calculated using 
the asymptotic expansion method. The theoretical prediction for the expected amplitude obtained using 
orthogonal basis and Floquet’s basis coincide. For the expected angular frequency, the Floquet’s basis gives 
a more accurate result, as a consequence of the partial decoupling between phase and amplitude. 

Finally figure Q] shows the stationary distribution for the probability density of the amplitude ( p(p,t) 
for t —> oo). The experimental probability to find the amplitude in the interval p + dp is evaluated as the 
fraction of simulation time spent in the interval. A Gaussian distribution with mean 1 + E[R] and variance 
var(i?) is also shown for comparison. 
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Figure 4: Stationary amplitude distribution p(p. t —> +oo), for the Stuart-Landau oscillator (1191) . A Gaussian distribution with 
mean // = ] + E[R] and variance a 2 = var(i?) is shown for comparison. Parameters are: a = 5, f) = 1, e = 0.15. 


6 . Conclusions 

A novel description for nonlinear oscillators subject to white Gaussian noise and described by Ito stochas¬ 
tic differential equations, is given. The dynamics is described in terms of a phase and an amplitude deviation 
variables. The phase function defines the projection of the stochastic orbit onto a reference trajectory, i.e. 
a limit cycle of the unperturbed system, while the amplitude represents the deviation from the reference 
orbit. 

The specific equations that are obtained depend upon the choice of a particular basis vectors. Two 
main basis have been considered: orthogonal basis and Floquet’s basis. Although the orthogonal basis 
is conceptually simpler, it has been shown that Floquet’s basis allows a partial decoupling between the 
amplitude and the phase equations. This permits to obtain more accurate estimation of the expected 
angular frequency, and to derive simplified phase reduced models. 

In the weak noise limit, the amplitude and phase equation can be conveniently solved with the help 
of asymptotic expansions. The method recast the equations in the form of a sequence of time dependent 
Ornstein-Uhlenbeck processes, that can be solved iteratively. In alternative, useful information concern¬ 
ing the expectation values of amplitude and angular frequency can be obtained using Ito calculus. Taking 
stochastic expectations and using the zero expectation property of Ito integrals a system of ordinary differ¬ 
ential equations for expectation values can be obtained. Due to the fact that Ornstein-Uhlenbeck processes 
are linear stochastic differential equations, the system is closed and can be solved analytically thus allowing 
a complete characterization of the stochastic process. 
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7 . Appendix A: Floquet theory basics 

Consider the n-dimensional homogeneous linear system of ODEs 

^ (78) 

where A[t ) is a nx n dimensional matrix with periodic entries. Let y ± (t ),..., y n (t) be n linearly independent 
solutions of (EH). Then 

• $(t) = [Ui(t),... ,y n (t)] is called a fundamental matrix. If $(0) = I n then <b(f) is called the principal 
fundamental matrix of the state transition matrix. 

• Any solution of (l78l) can be written in the form y(t) = <E>(t)y(0). 
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• The fundamental matrix is not unique, but they are all similar, that is, if $(f) and <F(t) are fundamental 
matrices, then a constant matrix C exists such that <h(f) = $>(t)C. 

• Let A(t) = A(t+T ), then $(f+T) is also a fundamental matrix, and $(i+T) = $(f).D. The eigenvalues 
of the constant, non singular matrix D, pi, ..., p n are called Floquet’s characteristic multipliers. They 
are related to the Floquet’s characteristic exponents, iq,..., v n by the formula pi = e UiT . 

Theorem 4. If y is a characteristic multiplier of system CHI)- then a non trivial solution y s (t) exists 
satisfying 

Vsit + T) =yy s (t) ( 79 ) 

Viceversa, if a nontrivial solution y s (t) of (1751) satisfies (1751) . then y is a characteristic multiplier and y s (0) 
is the corresponding eigenvector. 

Theorem 5 (Floquet 1883). The fundamental matrix 4>(f) of system (1781) can be written in the form 

$(t) = P(t)e Dt S 0 (80) 

where D = diag[v \,..., v n ], P(t) is a T-periodic regular matrix and So = P _ 1 (0) 

Theorem 6 . //4>(t) is a fundamental matrix of (1781) . then ’F(f) = [$“ 1 (t)] T is a fundamental matrix for 
the adjoint problem 

^± = -A T (t)z(t) (81) 

Therefore, denoting by p^t) the i th column of P(t) and by Si(t) the i th row of P _1 (t), if y s (t) = Pi(t) e Vit 
is a solution of (1751) . then z s (t) = Si(t) e~ Vit is a solution of (1511) . 


Appendix B: some details on Ito calculus 


Stochastic processes are nowhere differentiable. As a consequence the SDEs m should be interpreted as 
a shorthanded notation for the integral equation 

a(X s )ds+ I B(X S ) dW s (82) 

Jo 

Depending on the how the second integral on the right hand side is defined, different interpretations are 
possible. The most popular interpretations are Stratonovich and Ito. 

• Stratonovich integral. We adopt the standard notation f(X t ) o dW t to denote Stratonovich inter¬ 
pretation. In this view the functions are evaluated at the middle point of each time interval, that 
is 

rT n / x -L x \ 

l f(Xt) o dW t = ms -lirn^ ^ tj 2 (W tj - W tj _ x ) (83) 

where ms-lim denotes the mean square limit. 



• Ito stochastic integral. In Ito interpretation the functions are evaluated at the beginning of each 
interval 


f(X t )dW t = ms -lhnC*b-0 (W ti - W^) 


(84) 


3 = 1 


Ito lemma states that 


dt 2 = dt dWi = dWi dt = 0, dWi dWj = Sijdt 


(85) 
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• Ito SDEs do not follow the traditional calculus rule for change of variables. Ito formula must be 
used instead. Let X t be an Ito process solution of the SDE JTJ, and let g : R" H R" be a twice 
differentiable function, then g(X t ) is again an Ito process and its components are given by 


dg k = Y / J tdX i 


2—1 


Ox, 


1 n 


d 2 gk 


2 ^ dxidxi 


dXidX , 


( 86 ) 




Since in Ito interpretations stochastic processes and noise increments are independent, Ito integrals 
have zero expectation value 


E 



f(t)dW t 


= 0 


(87) 
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